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Abstract 


The second-order azimuthal anisotropy Fourier harmonics, Vj_, are obtained in pPb 
and PbPb collisions over a wide pseudorapidity (?/) range based on correlations 
among six or more charged particles. The pPb data, corresponding to an integrated 
luminosity of 35 nb^^, were collected during the 2013 LHC pPb run at a nucleon- 
nucleon center-of-mass energy of 5.02 TeV by the CMS experiment. A sample of semi¬ 
peripheral PbPb collision data at = 2.76 TeV, corresponding to an integrated 

luminosity of 2.5 and covering a similar range of particle multiplicities as the 
pPb data, is also analyzed for comparison. The six- and eight-particle cumulant and 
the Lee-Yang zeros methods are used to extract the V 2 coefficients, extending previ¬ 
ous studies of two- and four-particle correlations. For both the pPb and PbPb sys¬ 
tems, the V 2 values obtained with correlations among more than four particles are 
consistent with previously published four-particle results. These data support the in¬ 
terpretation of a collective origin for the previously observed long-range (large Arj) 
correlations in both systems. The ratios of V 2 values corresponding to correlations 
including different numbers of particles are compared to theoretical predictions that 
assume a hydrodynamic behavior of a pPb system dominated by fluctuations in the 
positions of participant nucleons. These results provide new insights into the multi¬ 
particle dynamics of collision systems with a very small overlapping region. 
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Measurements at the CERN LHC have led to the discovery of two-particle azimuthal corre¬ 
lation structures at large relative pseudorapidity (long-range) in proton-proton (pp) Jl] and 
proton-lead (pPb) EJEI collisions. Similar long-range structure has also been observed for 
= 200 GeV deuteron-gold (dAu) collisions at RHIC ||6HZl- The results extend previous 
studies of relativistic heavy ion collisions, such as for the copper-copper [(8|, gold-gold 15]- 
[l2l , and lead-lead (PbPb) ITOtiTSl systems, where similar long-range, two-particle correlations 
at small relative azimuthal angle \A(p\ « 0 were first observed. A fundamental question is 
whether the observed behavior results from correlations exclusively between particle pairs, or 
if it is a multi-particle, collective effect. It has been suggested that the hydrodynamic collective 
flow of a strongly interacting and expanding medium IIT9l - l2ll is responsible for these long- 
range correlations in central and mid-central heavy ion collisions. The origin of the observed 
long-range correlations in collision systems with a small overlapping region, such as for pp and 
pPb collisions, is not clear since for these systems the formation of an extended hot medium 
is not necessarily expected. Various theoretical models have been proposed to interpret the 
PP 12211231 arid pPb results, including initial-state gluon saturation without any final state in¬ 
teractions Il24l [25l and, similar to what is thought to occur in heavier systems, hydrodynamic 
behavior that develops in a conjectured high-density medium Il26] - l28l . These models have 
been successful in describing different aspects of the previous experimental results. 

To further investigate the multi-particle nature of the observed long-range correlation phe¬ 
nomena, in this Letter we present measurements of correlations among six or more charged 
particles for pPb collisions at a center-of-mass energy per nucleon pair of = 5.02 TeV. 

The azimuthal dependence of particle production is typically characterized by an expansion 
in Lourier harmonics (Vn) lEH. In hydrod}mamic models, the second (^ 2 ) and third (v^) har¬ 
monics, called "elliptic" and "triangular" flow 11301 , respectively, directly reflect the response 
to the initial collision geometry and fluctuations ||3T] - I33I , providing insight into the fundamen¬ 
tal transport properties of the medium. Lirst attempts to establish the multi-particle nature of 
the correlations observed in pPb collisions were presented in Refs. Il34ll35l by directly measur¬ 
ing four-particle azimuthal correlations, where the elliptic flow signal was obtained using the 
four-particle cumulant method 1361 . However, four-particle correlations can still be affected by 
contributions from non-collective effects such as fragmentation of back-to-back jets. By extend¬ 
ing the studies to six- and eight-particle cumulants II36II and by also obtaining results using the 
Lee-Yang zeros (LYZ) method, which involves correlations among all detected particles Il37ll38l . 
it is possible to further explore the collective nature of the correlations. High-statistics data ob¬ 
tained by the CMS experiment during the 2013 pPb run at the LHC are used. With a sample 
of very high final state multiplicity pPb collisions, the correlation data have been studied in 
a regime that is comparable to the charged particle multiplicity of the 50% most peripheral 
(semi-peripheral) PbPb collisions at = 2.76 TeV. 

The CMS detector comprises a number of subsystems [33 • The results in this Letter are mainly 
based on the silicon tracker information. The silicon tracker, located in the 3.8 T field of a su¬ 
perconducting solenoid, consists of 1440 silicon pixel and 15148 silicon strip detector modules. 
The silicon tracker measures charged particles within the pseudorapidity range |?/| < 2.5, and 
provides an impact parameter resolution of ~15 and a transverse momentum (px) reso¬ 
lution better than 1.5% at px ~ lOOGeV/c. The electromagnetic (ECAL) and hadron (HCAL) 
calorimeters are also located inside the solenoid and cover the pseudorapidity range |p| < 3.0. 
The HCAL barrel and endcaps are sampling calorimeters composed of brass and scintillator 
plates. The ECAL consists of lead tungstate crystals arranged in a quasiprojective geometry. 
Iron/quartz-fiber Cerenkov hadron forward (HE) calorimeters cover the range 2.9 < IpI < 5.2 
on either side of the interaction region. These HE calorimeters are azimuthally subdivided into 
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20° modular wedges and further segmented to form 0.175 x 0.175 rad {Arj x A(p) "towers". The 
detailed Monte Carlo (MC) simulation of the CMS detector response is based on Geant4 BoI . 

The analysis is performed using data recorded by CMS during the LHC pPb run in 2013. The 
data set corresponds to an integrated luminosity of 35nb~^. The beam energies were 4TeV 
for protons and 1.58 TeV per nucleon for lead nuclei, resulting in = 5.02 TeV. The beam 
directions were reversed during the run allowing a check of one potential source of systematic 
uncertainties. As a result of the energy difference between the colliding beams, the nucleon- 
nucleon center-of-mass in the pPb collisions is not at rest with respect to the laboratory frame. 
Massless particles emitted at rjcm = 0 in the nucleon-nucleon center-of-mass frame will be de¬ 
tected at rj = —0.465 (clockwise proton beam) or 0.465 (counterclockwise proton beam) in the 
laboratory frame. A sample of = 2.76 TeV PbPb data collected during the 2011 LHC 

heavy-ion run, corresponding to an integrated luminosity of 2.3 /rb^^, is also analyzed for com¬ 
parison purposes. The triggers and event selection, as well as track reconstruction and selection 
are summarized below and are identical to those used in Ref. |[35^ . 

Minimum bias (MB) pPb events were triggered by requiring at least one track with pj > 
0.4GeV/c to be found in the pixel tracker for a pPb bunch crossing. Only a small fraction 
(~10~^) of all MB triggered events were recorded (i.e., the trigger was "prescaled") because 
of hardware limits on the data acquisition rate. In order to select high-multiplicity pPb colli¬ 
sions, a dedicated high-multiplicity trigger was implemented using the CMS level-1 (LI) and 
high-level trigger (HLT) systems. At LI, three triggers requiring the total transverse energy 
summed over ECAL and HCAL to be greater than 20, 40, and 60 GeV were used since these 
cuts selected roughly the same events as the three HLT multiplicity selections discussed below. 
Online track reconstruction for the HLT was based on the three layers of pixel detectors, and 
required a track origin within a cylindrical region of length 30 cm along the beam and radius 
0.2 cm perpendicular to the beam around the nominal interaction point. For each event, the 
vertex reconstructed with the highest number of pixel tracks was selected. The number of pixel 
tracks with \r]\ < 2.4, pj > 0.4 GeV/c, and a distance of closest approach to this vertex 

of 0.4 cm or less, was determined for each event. Several high-multiplicity ranges were defined 
with prescale factors that were progressively reduced until, for the highest multiplicity events, 
no prescaling was applied. 

In the offline analysis, hadronic collisions are selected by requiring a coincidence of at least one 
HE calorimeter tower containing more than 3 GeV of total energy in each of the HE detectors. 
Only towers within 3 < |?/| < 5 are used to avoid the edges of the HE acceptance. Events are 
also required to contain at least one reconstructed primary vertex within 15 cm of the nominal 
interaction point along the beam axis and within 0.15 cm transverse to the beam trajectory. At 
least two reconstructed tracks are required to be associated with the primary vertex. Beam re¬ 
lated background is suppressed by rejecting events for which less than 25% of all reconstructed 
tracks pass the track selection criteria of this analysis. The pPb instantaneous luminosity pro¬ 
vided by the EHG in the 2013 run resulted in an approximately 3% probability of at least one 
additional interaction occurring in the same bunch crossing. Following the procedure devel¬ 
oped in Ref. Il35l for rejecting such "pileup" events, a 99.8% purity of single-interaction events is 
achieved for the pPb collisions belonging to the highest multiplicity class studied in this Fetter. 
In pPb interactions simulated with the EPOS Il4l]l and HIJING II42I event generators, requiring at 
least one primary particle with total energy E > 3 GeV in each of the tj ranges — 5 < ?/ < — 3 
and 3 < p < Bis found to select 97-98% of the total inelastic hadronic cross section. 

The GMS "high-quality" tracks, described in Ref. Il43l are used in this analysis. Additionally, 
a reconstructed track is only considered as a candidate track from the primary vertex if the 
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significance of the separation along the beam axis (z) between the track and the best vertex, 
dz/cr{dz), and the significance of the track-vertex impact parameter measured transverse to the 
beam, dj/cr{dj), are each less than 3. The relative uncertainty in the transverse-momentum 
measurement, 0'{pj)/pj, is required to be less than 10%. To ensure high tracking efficiency 
and to reduce the rate of incorrectly reconstructed tracks, only tracks within \ t]\ < 2.4 and with 
0.3 < Pj < 3.0 GeV/c are used in the analysis. A different pj cutoff of 0.4GeV/c is used in the 
multiplicity determination because of constraints on the online processing time for the HLT. 

The entire pPb data set is divided into classes of reconstructed track multiplicity, The 

multiplicity classification in this analysis is identical to that used in Ref. Il35ll . where more de¬ 
tails are provided including a table relating to the fraction of MB triggered events. A 

subset of semi-peripheral PbPb data collected during the 2011 LHG heavy-ion run with an MB 
trigger is also reanalyzed in order to directly compare the pPb and PbPb systems at the same 
track multiplicity. This PbPb sample is reprocessed using the same event selection and track re¬ 
construction as for the present pPb analysis. A description of the 2011 PbPb data can be found 
in Ref. [I44i|. 

Extending the previous two- and four-particle azimuthal correlation measurements of Ref. Il35l . 
six- and eight-particle azimuthal correlations II36II are evaluated in this analysis as: 


= ^(^Qin((pl+(p2+4’3+‘l>i-<p5-<p6-<p7- 


( 1 ) 


Here (pi {i = 1, ..., 8) are the azimuthal angles of one unique combination of multiple particles in 
an event, n is the harmonic number, and ((•••)) represents the average over all combinations 
from all events within a given multiplicity range. The corresponding cumulants, c„{6} and 
c„{8}, are calculated as follows: 

Cn{6}=«6»-9.{{4)>{{2))+12.{{2))=, 

^,{8} =({8)> - 16 . ({6))((2)) - 18 . ((4))^ + 144 - {{4))({2))2 - 144{{2))‘, 


using the Q-cumulant method as formulated in Ref. Il36l , where ((2)) and ((4)) are defined 
similarly as in Eq. Q. The Eourier harmonics Vn that characterize the global azimuthal behavior 
are related to the multi-particle correlations 1451 using 


^ 46 } = 

= y-^c„{8}. 


(3) 


To account for detector effects, such as the tracking efficiency, the Q-cumulant method was 
extended in Ref. Il45l to allow for particles having different weights. Each reconstructed track is 
weighted by a correction factor to account for the reconstruction efficiency, detector acceptance, 
and fraction of misreconstructed tracks. This factor is derived as a function of pj and rj, as 
described in Refs. IT^ [HI, based on MG simulations. The combined geometrical acceptance 
and efficiency for track reconstruction exceeds 60% for pj 0.3GeV/c and \t]\ < 2.4. The 
efficiency is greater than 90% in the \rj\ < 1 region for pj > 0.6 GeV/c. Eor the entire multiplicity 
range (up to 350) studied in this Letter, no dependence of the tracking efficiency on 

multiplicity is found and the rate of mis-reconstructed tracks remains at the 1-2% level. The 
software package provided by Ref. 11451 is used to implement the weights of individual tracks 
in the cumulant calculations. 
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The LYZ method |l3Zl[38l allows a direct study of the large-order behavior by using the asymp¬ 
totic form of the cumulant expansion to relate locations of the zeros of a generating function 
to the azimuthal correlations. This method has been employed in previous CMS PbPb analy¬ 
ses OZlilSl. For each multiplicity bin, the V 2 harmonic averaged over 0.3 < pj < 3.0GeV/c is 
found using an integral generating function IHTI . Similar to the cumulant methods, a weight 
for each track is implemented to account for detector-related effects. In both methods, the 
statistical uncertainties are evaluated from data by dividing the data set into 20 subsets with 
roughly equal numbers of events and evaluating the standard deviation of the resulting distri¬ 
butions of the cumulant or U 2 {LYZ} values. In the case of low multiplicity or small flow signal, 
the LYZ method may overestimate the true collective flow. This effect was studied using MC 
pseudo-experiments for the event multiplicities covered in this analysis and a small correction 
is applied to the data. The correction is less than 3% in the lowest multiplicity bin and becomes 
much smaller in higher-multiplicity bins. This correction is also included in the quoted LYZ 
systematic uncertainties. 

Systematic uncertainties are estimated by varying the track quality requirements, by compar¬ 
ing the results using efficiency correction tables from different MC event generators, and by 
exploring the sensitivity of the results to the vertex position and to the bin width. For 

the pPb data, potential HLT bias and pileup effects are also studied by requiring the presence of 
only a single reconstructed vertex. No evident or beam direction dependent systematic 

effects are observed. For pPb collisions, a 5% systematic uncertainty is obtained for V 2 { 6 } and 
a 6 % uncertainty is found for both ! 72 { 8 } and U 2 {LYZ}. The corresponding uncertainties for 
PbPb collisions are 2% for ? 72 { 6 } and r’ 2 { 8 }, and 4% for ? 72 {LYZ}. 

In Fig.[^ the six- and eight-particle cumulants, 02 ( 6 } and C 2 { 8 }, for particle pj of 0.3-3.0GeV/c 
in 2.76 TeV PbPb and 5.02 TeV pPb collisions are shown as a function of event multiplicity. 
The cumulants shown are required to be at least two standard deviations away from their 
physics boundaries ic 2 { 6 }/ 0 'c 2 { 6 } > ^ 2 { 8 }/o'c 2 { 8 } < — 2 ), so that the statistical uncertainties 
can be propagated as Gaussian fluctuations [^I. Non-zero multi-particle correlation signals are 
observed in both PbPb and pPb collisions. The pPb data exhibit larger statistical uncertainties 
than the PbPb results, mainly because of the smaller magnitudes of the correlation signals. 
Because of the limited sample size, the 02 ( 6 } and 02 ( 8 } values in pPb collisions are derived for 
a smaller range in N°®™. 

The second-order anisotropy Fourier harmonics, V 2 , averaged over the px range of 0.3-3.0 GeV/c, 
are shown in Fig. [phased on six- and eight-particle cumulants (Eq. for 2.76 TeV PbPb (left) 
and 5.02 TeV pPb (right) collisions, as a function of event multiplicity. The open symbols are V 2 
results extracted by GMS using two- and four-particle correlations 1351 . The V 2 values derived 
using the LYZ method involving correlations among all particles are also shown. For each mul¬ 
tiplicity bin, the values of r' 2 { 4 }, V 2 { 6 }, i^ 2 { 8 }, and i? 2 {LYZ} for pPb collisions are found to be 
in agreement within 10%. For part of the multiplicity range, the values for i> 2 { 4 } are larger than 
the others by a statistically significant amount, although still within 10%. The corresponding 
PbPb values are consistently higher than for pPb collisions, but within the PbPb system are 
found to be in agreement within 2 % for most multiplicity ranges and within 10 % for all mul¬ 
tiplicities. This supports the collective nature of the observed correlations, i.e., involving all 
particles from each system, and is inconsistent with a jet-related origin involving correlations 
among only a few particles. The V 2 data from two-particle correlations are consistently above 
the multi-particle correlation data. This behavior can be understood in hydrodynamic mod¬ 
els, where event-by-event participant geometry fluctuations of the V 2 coefficient are expected 
to affect the two- and multi-particle cumulants differently Il48ll4^ . Note that, to minimize jet- 
related non-flow effects, the 1 ^ 2 { 2 } values are obtained with an rj gap of 2 units between the 
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Figure 1: The cumulant 02 ( 6 } and — 02 ( 8 } results as a function of for PbPb and pPb 

reactions. Error bars and shaded areas denote statistical and systematic uncertainties, respec¬ 
tively 

two particles. Possible residual non-flow effects resulting from back-to-back jet correlations are 
estimated using very low multiplicity events in Ref. I35l . Based on this analysis, such non-flow 
effects are expected to make a negligible contribution to V 2 { 2 } in very high multiplicity events. 
In PbPb collisions, the V 2 values from all methods show an increase with multiplicity, while lit¬ 
tle multiplicity dependence is seen for the pPb data. This difference might reflect the presence 
of a lenticular overlap geometry in PbPb collisions, which is not expected in pPb collisions, that 
gives rise to a large (and varying) initial elliptic asymmetry in the PbPb system. 

The effect of fluctuation-driven initial-state eccentricities on multi-particle cumulants has re¬ 
cently been explored in the context of hydrodynamic behavior of the resulting medium ||50l|5Tl. 
For fluctuation-driven initial-state conditions, ratios of V 2 values derived from various orders 
of multi-particle cumulants are predicted to follow a universal behavior IISOl . In Fig.l^ ratios of 
r' 2 { 6 }/f 2 { 4 } (top) and U 2 { 8 }/r' 2 { 6 } (bottom) are calculated and plotted against U 2 ^}/i’ 2 { 2 } 
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Figure 2: The V 2 values as a function of Open data points are published two- and four- 

particle V 2 results Il35l . Solid data points are V 2 results obtained from six- and eight-particle 
cumulants, and LYZ methods, averaged over the particle pj range of 0.3-3.OGeV/c, in PbPb at 
= 2.76 TeV (left) and pPb at = 5.02 TeV (right). Statistical and systematic uncertain¬ 
ties are indicated by the error bars and shaded regions, respectively. 


in pPb collisions at = 5.02 TeV. The ^ 2 ( 2 } and U 2 { 4 } data are taken from previously pub¬ 
lished CMS results ||35|| . The solid curves correspond to theoretical predictions for both large 
and small systems based on hydrodynamics and the assumption that the initial-state geome¬ 
try is purely driven by fluctuations [501. The ratios from PbPb collisions are also shown for 
comparison. Note that the geometry of very central PbPb collisions might be dominated by 
fluctuations, but for these semi-peripheral PbPb collisions the lenticular shape of the overlap 
region should also strongly contribute to the U 2 values. The CMS pPb data are consistent with 
the predictions within statistical and systematic uncertainties. The systematic uncertainties in 
the ratios presented in Fig. are estimated to be 2.4% for U 2 { 4 }/u 2 { 2 } for both pPb and PbPb 
collisions, 1% for r' 2 { 6 }/r' 2 { 4 } in pPb and PbPb collisions, and 3.6% and 1% for V 2 { 8 }/V 2 { 6 } 
in pPb and PbPb collisions, respectively. Since they are all derived from the same data, the 
systematic uncertainties for the different cumulant orders are highly correlated and therefore 
partially cancel in the ratios. 

Recently, other theoretical models based on quantum chromodynamics, and not involving hy¬ 
drodynamics, have also been suggested to explain the observed multi-particle correlations in 
pPb collisions Il52ll53ll . Unlike the descriptions based on hydrod}mamic behavior, these models 
do not require significant final-state interactions among quarks and gluons. They suggest sim¬ 
ilar values for U 2 { 4 }, U 2 { 6 }, U 2 { 8 }, and U 2 {LYZ}, without yet, however, providing quantitative 
predictions. 

In summary, multi-particle azimuthal correlations among six, eight, and all particles have been 
measured in pPb collisions at y^s^ = 5.02 TeV by the CMS experiment. The new measure¬ 
ments extend previous CMS two- and four-particle correlation analyses of pPb collisions and 
strongly constrain possible explanations for the observed correlations. A direct comparison of 
the correlation data for pPb and PbPb collisions is presented as a function of particle multi- 
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Figure 3: Cumulant ratios r’2{6}/u2{4} (top) and V 2 { 8 }/v 2 { 6 } (bottom) as a function of 
1^2 {4}/ 1 ^ 2 {2} inpPb collisions at = 5.02 TeV and PbPb collisions at = 2.76 TeV. Prror 
bars and shaded areas denote statistical and systematic uncertainties, respectively. The solid 
curves show the expected behavior based on a hydrod}mamics motivated study of the role of 
initial-state fluctuations |[50l . 

plicity. Averaging over the particle pj range of 0.3-3.0 GeV/c, multi-particle correlation signals 
are observed in both pPb and PbPb collisions. The second-order azimuthal anisotropy Fourier 
harmonic, V 2 , is extracted using six- and eight-particle cumulants and using the LYZ method 
which involves all particles. The V 2 values obtained using correlation methods including four 
or more particles are consistent within ±2% for the PbPb system, and within ±10% for the pPb 
system. This measurement supports the collective nature of the observed correlations. The ra¬ 
tios of V 2 values obtained using different numbers of particles are found to be consistent with 
hydrodynamic model calculations for pPb collisions. 
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